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ABSTRACT 

An approach, based on a global averaging procedure, is presented 
for estimating the power spectrum of a second order stationary 
zero- mean ergodic stochastic process from a finite length record. 
This estimate is derived by smoothing, with a cubic smoothing 
spline, the naive estimate of the spectrum obtained by applying FFT 
techniques to the raw data. By means of digital computer simulated 
results, a comparison is made between the features of the present 
approach and those of more classical techniques of spectral estima- 
tion. 
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INTRODUCTION 


Let [X « 3 <t<“} denote a second order stationary zero- 

average ergodic real-valued stochastic process. If the process is Gaussian, 

then 4t ds-completely characterized by the autocovariance' function 

C(t)= E[X(t) X(t + r)] 
l T/2 

= lim ) x(t) x(t + T)dt, (1. 1) 

T - « 1 -T/2 

or, equivalently, by the Fourier transform of C(t) , referred to as the 
power spectrum, 


P(f) = J C(t) e" i2nfT dT 

— CO 


T/2 


= lim 

P ~t CO 


T 


x(t) 

/O 


~i2nft 


dt 


-T/2 


( 1 . 2 ) 


Above and henceforth, we denote by x(t) or x^. a realization of X(t) 


or X t . 


If it is desired to estimate P(f) on the basis of a record of 


finite length T^ , we may compute the apparent autocovariance function 


1 

C 00 ( T ) = T - t| J x ( t_ / 2 ) x ( t+ /o) 
n 1 -b 


dt 


(1. 3) 


where b = (T r - | t | / 2 ) 
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and take the Fourier transform of Cqq(t) = Cqq(t) for | t j<T m < T n , = 0 
otherwise, giving 

co 

? 00 (t) = I <1> 9)0 (T)e " i2TlfTdT • (1- 4) 

/v< 

Unfortunately, in general Pqq( 1) does not converge to P(f) , in 
the -norm, as the length of the time series becomes infinite, i. e. , 
P 00 (f) is not a consistent estimator for P(f). 

The Tukey family of filters^ remedies this difficulty by an 
averaging process 

ITT 

C a ( t ) = (a -Kl - a) cos C 00 (t) 


(1. 5) 


for t sT < T , 
II m n 

--- 0 if It! > T , 

I 1 m 


where 0 < a < 1 , giving 

p a< t > = J* 


/ _ \ “ i 2 i if t , 

C a ( t ) e dT 


(1. 6) 


For the case where the record { x(t): t e (0, T^)} is sliced at n intervals of 
length At, we first compute the mean lagged products 


1 


C =- 


q=n- r 

r n - r x q x q+rh (1. 7) 

for r = 0,l,2, . . . , m where m< n. We then apply a discrete finite cosine 


series transform to C Q , . . . , G m obtaining 


V A t 
r = 




m 
m-1 




C 0 + 2 I C q cos 


qrTT , „ 

~ — + C cos rn 
rn m 




q=l 


(1. 8) 


J 
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The Tukey filter then gives as estimates of the spectral density for the case 
where 

1 

a = 2 (hanning) 


U Q = . 5V q + . 5V 1 


u r - ■ 25VJ..J + . 5V r + . 25V r+J , lsrsm-1 


U m = - 5V m-l + - 5V m- 


Basically, then, the popular Tukey filter is a local averaging procedure. 

In this paper we present, as an alternative approach, an estimation proce- 
dure based on a global averaging approach. The results given here are of a 
preliminary nature and are now in the process of being studied in greater 
detail. 


2. DISCUSSION AND RESULTS 
Consider the discrete white noise time series 

. 

X t = e t ; t = 1,2, , . . , n , . (2. 1) 

where the e are independent Gaussian random variables with mean zero 
and variance unity. A naive estimate of the power spectrum is given by 

X x t e' i2TTft 2 (2. 2) 

t=l 

and may be quickly computed by Fast Fourier Transform techniques^. 

A plot of this estimate is given for a time series of length n = 512 in 
Figure 1. The inconsistency of this estimate is apparent. We shall attempt 
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to remedy this inconsistency by smoothing 


P N< f > 


Let us divide the frequency interval [ 0, . 5cps] into k-1 equal 
subintervals. We seek the solution P (f) which minimizes , V P e C 1 [0, .5] 

b O 

and piecewise continuous P " (See (3) ), 

o 


k -1 

z 

j =0 


r 


a. 

J 




J 


s v 2 (k-l) 


) -Pxi( 


J 


N v 2(k-1) 


') 


+ g 


J 


. 5 


(P s ’ (f)) df , 

V y < 2 ’ 3) 

where 5 is a nonnegative constant which can be selected so as to give 
appropriate weight to the roughness penalty term, and the a^ are positive 
constants summing to unity used for weighting the "miss distances" in the 
first term. In the ensuing discussion we shall take 


1 

a j = k '» for 3 = 0 , 1 , 2 , ...» k -1 . 

Schoenberg has shown ^ that the solution P is, in fact, a cubic 

spline with knots st the subinterval endpoints. The smoothing algorithm 

employed here is that of Greville^). 

In Figure 2, we show P (f) plots for k = 8 and n = 128 and 

s 

512 together with the true power spectrum. We note a dramatic change 
from the naive estimate and a good approximation to the actual spectral 
density. 


Next, let us consider a spectral estimation problem from Jenkins 

and Watts' Spectral Analysis and its Applications^’ PP* ^3 ^65). jj ie 

2 

spectral density P(f)/cr of the first order autoregressive process 

X 

X^. = - . + Sj. ; t = 1, 2, . . . , n , (2. 4) 
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where the e are independent Gaussian variates with zero mean and unit 
variance, is given by 

2 1. 68 1 
p (f)A? x = 1. 16 + . 8 cos 2nf ’ 0 * f ^ 2 ' 5) 

In Figure 3, we demonstrate estimates of the spectral density for series of 

lengths 128, 256 and 512 using hanning with m = 32 . In Figure 4, we show 

'P (f) for the same time series, 16 knots and g - 1 . The spline smoothing 
s 

approach appears to compare favorably with hanning. 

The effect of the number of knots is demonstrated in Figure 5 in 
the estimation of the spectral density of the autoregressive process in (2. 4) 
for a series of length 1024 using g = 1 and k = 8, 16, 32. In Figure 6, we 
demonstrate estimates using the same series with 16 knots and g = 1, .1, 
and .01. Clearly, the selection of g and k is an important consideration 
and should be based both on the length of the time series and prior feelings 
as to the shape of the spectral density. We have not yet had opportunity to 
investigate optimal (g, k) selection procedures but shall approach this 
problem in the future. 

As a final example from Jenkins and Watts ^68 272)^ we 

consider the second order autoregressive process 

X t = X t _ 1 - . 5X t _ 2 + e t ; t = 1, 2, . . . , n , (2. 6) 

where the e are independent Gaussian deviates with mean zero and unit 

variance. Tie spectral density here is given by 

2 . 834 1 

P(f)/ CT x " 2. 25-3 cos 2nf + cos 4irf ’ ^ ^ ^ ^2 


(2. 7) 
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In Figure 7, we show for a series of length 512, the hanning estimate with 
m = 32 and compare the result with the spline smoothing approach using 
g = .01 and 32 knots. 

3. CONCLUSIONS 

We have explored in a preliminary fashion the possibility of 
using spline smoothing techniques, in the frequency domain, for the esti- 
mation of power spectra. Although for the examples considered, this 
technique appears to compare favorably with one of the most popular 
smoothing techniques -hanning - we still need to study in detail the consis- 
tency properties of the proposed estimation procedure. We hope, in addi- 
tion, to study optimal methods for applying spline smoothing to cross - 
spectra. 
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Fig. 2: Spline smoothed spectral density 
estimates for white noise. 
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Fig. 3: Hanning smoothed spectral density- 
estimates for a first-order ar 
process = a X^ + e • 


Fig- 4: Spline smoothed spectral densi 
estimates for a first order ar 
process X = a x , + E 
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Fig. 5: Spline smoothed spectral density 

estimates for a first order ar 
process X = X ^ under 

different conditions on knots. 


Fig. 6: Spline smoothed spectral density 

estimates of a first order ar process 
under different conditions on the 
smoothing parameter g . 
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Fig. 7: Hanning and spline smoothed 

estimates of the spectral density 
of a second order ar process 
X = a. X , + a- ^ + £, 


